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In this paper we present an accurate numerical scheme for extracting inter-atomic exchange param¬ 
eters (Jij) of strongly correlated systems, based on first-principles full-potential electronic structure 
theory. The electronic structure is modelled with the help of a full-potential linear muffin-tin orbital 
(FP-LMTO) method. The effects of strong electron correlations are considered within the charge 
self-consistent density functional theory plus dynamical mean-field theory (DFT+DMFT). The ex¬ 
change parameters are then extracted using the magnetic force theorem, hence all the calculations 
are performed within a single computational framework. The method allows to investigate how the 
Jij -parameters are affected by dynamical electron correlations. In addition to describing the for¬ 
malism and details of the implementation, we also present magnetic properties of a few commonly 
discussed systems, characterised by different degrees of electron localisation. In bee Fe, treated as 
a moderately correlated metal, we found a minor renormalisation of the Jij interactions once the 
dynamical correlations are introduced. However, generally, if the magnetic coupling has several com¬ 
peting contributions from different orbitals, the redistribution of the spectral weight and changes in 
the exchange splitting of these states can lead to a dramatic modification of the total interaction 
parameter. In NiO we found that both static and dynamical mean-field results provide an adequate 
description of the exchange interactions, which is somewhat surprising given the fact that these 
two methods result in quite different electronic structures. By employing Hubbard-I approximation 
for the treatment of the 4/ states in hep Gd we reproduce the experimentally observed multiplet 
structure. The calculated exchange parameters result to be rather close to the ones obtained by 
treating the 4/ electrons as non-interacting core states. 


INTRODUCTION 

The interatomic exchange interaction (Jij), together 
with the magnetic moment, are the key quantities for a 
microscopic description and understanding of magnetism 
in real materials. Being of a purely quantum origin, they 
define most of the macroscopic properties of magnetic 
materials, in particular their Curie temperature (T c ), the 
magnon dispersion and the magnetic response to exter¬ 
nal stimuli such as an applied magnetic field. Having a 
sufficiently large T c , preferably above room temperature, 
is the key ingredient for any magnetic material that is to 
be used in the new generation electronic devices (see e.g. 
Ref. 1), or for applications in more traditional fields of 
magnetism such as new permanent magnetic materials. 
Therefore, an ability to predict from first principles elec¬ 
tronic structure theory the exchange couplings in various 
materials is an essential step towards achieving important 
goals in the design of materials for emerging technologies. 

From a theoretical viewpoint, the basic microscopic 
mechanisms which give rise to the exchange couplings 
are quite well-known. 2 However, when it comes to real 
systems, predicting the sign and the magnitude of the 
-parameters is a challenging theoretical problem. In 
fact the electronic structure of real materials is complex 
and results into a competition between various exchange 


mechanisms. Which one happens to dominate is rather 
difficult to judge a priori. 

In 1987, a general scheme was proposed for the extrac¬ 
tion of the exchange integrals from electronic structure 
calculations. 3 (The basic idea was published few years 
earlier in Ref. 4) Within this approach the energy of the 
electronic Hamiltonian is mapped onto a classical Heisen¬ 
berg model of the form: 

H = — JijOi • ej, (1) 

iAj 

where denotes the unit vector along the magnetic mo¬ 
ment at the site i. The Jij -parameters between sites i 
and j are defined as the response to infinitesimally small 
rotations of the corresponding atomic spins, i.e. in a 
linear-response manner. The approach was successfully 
applied to a variety of different materials and has be¬ 
come the primary tool for studying inter-site magnetic 
interactions from first principles. 5 ~ 10 Note that the ex¬ 
change parameters defined in this way are, in general, 
dependent on initial magnetic configuration. In itinerant 
systems such a dependence may be quite strong, which 
implies the non-Heisenbergian character of the magnetic 
interactions. 11 

In 2000 it was shown that in the presence of dynam¬ 
ical correlations the formalism is still valid if the self- 
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energy E is local. 12 This is exactly the case of the Dy¬ 
namical Mean-Field Theory (DMFT) 13 , which, combined 
with the density functional theory (DFT), is the state-of- 
the-art method for ab initio modelling of the electronic 
structure of correlated materials. 14-16 In 2006 X. Wan 
et al used the same formalism to study the prototypi¬ 
cal materials with strong electronic correlations, i.e. the 
transition-metal monoxides. 17 Despite Wan and cowork¬ 
ers employed several approximations to make the calcula¬ 
tions feasible, they managed to successfully compute the 
exchange parameters, which were then used to simulate 
the magnon spectra. The authors also reported that 
the best agreement between theory and experiment was 
achieved by employing a DMFT-type technique. 

In this paper we describe our implementation of the 
formalism of Ref. 12 in a full-potential electronic struc¬ 
ture method, in which strong electron correlations are 
treated with DMFT 18-23 . We apply this implementa¬ 
tion to several materials with different strength of cor¬ 
relations. After this Introduction, the present article is 
organised in the following way. In section I the general 
computational scheme is presented. Section II is ded¬ 
icated to the discussion about how the choice of local 
projections affects the computed exchange parameters. 
Sections III, IV and V report the results obtained for 
bee Fe, NiO, and hep Gd. The last section will highlight 
the main conclusions of our study, and perspectives of 
our future research. 


I. THEORETICAL METHODS 


We model the electronic structure with the help of the 
full-potential linear muffin-tin orbitals (FP-LMTO) code 
“RSPt”. 18,19 In order to treat strongly correlated materi¬ 
als, where the conventional DFT fails, we use a combina¬ 
tion of DFT and DMFT (DFT+DMFT), as implemented 
in a FP-LMTO electronic structure method. Details of 
this implementation have been presented elsewhere 20-23 , 
and will not be repeated here with an exception of some 
specific details that will be needed for our discussion. In 
DFT+DMFT one first selects a subset of electrons that 
are not described properly by standard local or semi-local 
exchange-correlation functionals. Then the Kohn-Sham 
Hamiltonian of the DFT problem is corrected with an 
explicit two-particle term describing the local Coulomb 
repulsion U among these localized electrons minus a dou¬ 
ble counting term. The latter assures that the contri¬ 
butions of added Coulomb term are removed from the 
Kohn-Sham Hamiltonian, where they are described less 
accurately (due to the failure of the exchange-correlation 
functional). The obtained Hamiltonian is an effective 
Hubbard model and can be solved by DMFT. Once the 
DMFT cycle is converged, a new electron density is ob¬ 
tained, which leads to a new Kohn-Sham Hamiltonian. 
All simulations performed in this study are fully con¬ 
verged with respect to the electron density. 20 

The effective impurity problem, which is at the core of 


the DMFT method, can be solved by means of various 
techniques. The results presented in Section III, for bee 
Fe, were obtained through the Spin-polarized T-matrix 
plus Fluctuation Exchange (SPTF) solver 24 , which is 
based on perturbation theory. The solver enforces the 
system to have the Fermi liquid properties and has to 
be used with care, i.e. only when U is smaller than the 
bandwidth W . This is usually true for the transition 
metals. The results presented in Section IV for NiO, 
instead, were obtained by means of the Exact Diagonal- 
isation (ED) solver. This is one of the most accurate 
techniques and can be applied for an arbitrary correla¬ 
tion strength. Recently ED has been applied to study 
the spectral properties of the transition metals monox¬ 
ides, and lead to a great agreement between theory and 
experiment. Finally the results presented in Section VII, 
for hep Gd, were obtained in the Hubbard I approxima¬ 
tion (HIA), 15 where the hybridisation between the cor¬ 
related orbitals and the rest is neglected, hence the im¬ 
purity electrons can only redistribute themselves within 
a given /-shell. The HIA approximation is usually em¬ 
ployed for compounds with extreme localisation of cer¬ 
tain valence electrons. For example, it was applied to a 
series of rare-earth-based compounds and showed a good 
agreement with to the measured photoemission spectra 
and ground state properties. 25-28 

Once the electronic structure problem is solved, one 
can extract the couplings. To this aim we should de¬ 
fine a set of local orbitals |i, £, a) for the atomic site i. In 
a crystal i indicates both the lattice vector R of the unit 
cell and the basis vector r within the unit cell, so that it 
points uniquely to a given ion. The exact shape of the 
local orbitals does not need to be specified now, and we 
will postpone this topic to the next Section. For now we 
keep the formulation as general as possible, and we use 
the global index £ to refer to a set of quantum numbers 
labelling the basis functions, a is a spin index, and can 
have values With this in mind, we can formu¬ 

late a generalized expression for the inter-site exchange 
parameters in DFT+DMFT: 


Jij — a ^ ^ Tr (icj n )C —(icj n ) Aj (icj n )C^-(icj n ) 


( 2 ) 


where the trace is intended over the orbital degrees of 
freedom £, T is the temperature, and uj n = 27rT(2n + l) is 
the n-th fermionic Matsubara frequency. G is the inter¬ 
site Green’s function between sites i and j and projected 
over a given spin a. We use only one spin index, since 
in this whole article we assume that contributions due to 
the spin-orbit coupling are negligible. The term A^ gives 
the exchange splitting at the site i: 

A i(iu n ) = H^l + £- H^i - £+w„), (3) 

where we introduced the spin- and site-projected Kohn- 
Sham Hamiltonian H ks and local self-energy E. The lat¬ 
ter comes from the solution of the DMFT equations. In 
absence of self-energy we recover the standard DFT ex¬ 
pression 3 , with a static exchange splitting. In DMFT, 12 
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instead, due to the fact that the self-energy is dynamical, 
the exchange field becomes dynamical too. One has to re¬ 
alise that E, in spite of being local, also affects the inter¬ 
site Green’s function. This is so, because the Green’s 
function is the resolvent of the entire Hamiltonian H : 


G{z) 


1 

z-H 


(4) 


Hence the inversion operation mixes all the matrix ele¬ 
ments: site on- and off-diagonal ones. From a numerical 
viewpoint, the computation of Eq. (2) is not demanding 
at all, even if the sum runs over an infinite number of 
Matsubara frequencies. Since each Green’s function at 
high frequencies behaves as ^ 1 /iuj n , the quantity under 
summation decays at least as fast as ^ 1 /(cj 2 ). In prac¬ 
tice the amount of Matsubara points used for evaluation 
of the self-energy in a standard DMFT cycle (the exact 
number depends on the temperature T) is more than suf¬ 
ficient for a reliable evaluation of the -parameters. 


II. LOCAL ORBITALS AND EXCHANGE 
PARAMETERS 

In this section we provide a detailed description of 
the construction of the localised basis {|z,£,cr)} for the 
evaluation of the exchange couplings. Given that also 
the DFT+DMFT method requires the definition of local 
orbitals, we will first present our formalism referred to 
this method. Later on, we will contextualize the dis¬ 
cussion to the exchange parameters. The LMTO ba¬ 
sis functions 29 that are used to solve the DFT problem 
are {|i?, r, %, <r)}, where R and r refer respectively to 
the lattice site and the basis vector within the unit cell, 
while a is again the spin. From {|R, r, X,cr)} one can 
construct the Bloch sums {|fc, r, y, a)} used to solve the 
DFT eigenvalue problem and subsequently for the DMFT 
one-electron Green’s function. In FP-LMTO the index x 
stands for {n, /, m, ft}, where n is the principal quantum 
number, l the angular quantum number, m the magnetic 
quantum number, and k the value of the kinetic energy 
tail, which defines the asymptotic behaviour of the func¬ 
tions at large distances. It is possible to use a minimal 
basis set, using only a single set of n,/,m and n. How¬ 
ever, for a better converged basis, it is useful to adopt a 
basis with several k values. This ensures the flexibility 
of the basis set to describe various shapes of the electron 
density (see Ref. 19 for more details). 

We consider two types of local orbitals {|i, £, cr)}, which 
we label as “ORT” and “MTH”. These orbitals are 
schematically represented in Fig. 1. The localized or¬ 
bitals “ORT” are constructed directly from the LMTO 
basis functions {|R, r, x,cr)} (also depicted in Fig. 1) 
through the following Lowdin orthogonalisation: 

MMobx = J2 eikR ' 

k,X 


(a) Original LMTO basis 



FIG. 1. Schematic view of the various local orbitals used in 
RSPt: (a) the non-orthonormal LMTO basis set used in DFT; 
(b) the Lowdin orthonormalised LMTO “ORT” used as first 
projection scheme; (c) the muffin-tin head projected orbitals 
“MTH” orbitals used as second projection scheme. 


where S £ ^ is the overlap matrix between the LMTO 
Bloch sums and the local orbitals in the unit-cell at the 
origin Rq. Hence, the (|i,£,cr)} are obtained from the 
original basis by performing a /c-point-wise orthonormal¬ 
isation. In fact, such transformation mixes all the wave- 
functions, including the ones centered at different atomic 
sites. Such mixing is anyway small, and the whole con¬ 
struction is analogous to the construction of the Wannier 
functions. 30 In practice, due to the long decaying tails of 
the LMTO, these orbitals are not very localized, which 
is more physical from the electronic point of view. More¬ 
over, there is a one-to-one correspondence between the 
index £ and the index y. This may be a problem if one 
has DFT+DMFT simulations in mind, since an atomic- 
like basis is required for a proper parametrization of the 
Coulomb interaction. For the ORT basis this condition is 
satisfied only when using a single k value. However, mul¬ 
tiple k values are often needed for a proper description 
of the interstitial region between atomic spheres. 

The localized orbitals “MTH”, on the other hand, re¬ 
sult from a projection onto the so-called muffin-tin head 
of the LMTO. The MTH functions are written as: 


( r I h £> (7 )mth 



\ri\ < Rmt^^ 
otherwise, 


where r* = r — Ri and <f> reads as the solution of the 
radial Schrodinger equation inside the MT sphere at the 
linearization energy 22 . The MTH orbitals are extremely 
localised, since they are basically atomic functions. They 
have a pure /,m character and do not overlap by con- 
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struction. For DFT+DMFT the MTH set is particularly 
convenient, as in this case £ = {n,/,ra}, as clear from 
Eq. (6). There is no need to be in a one-to-one rela¬ 
tion with x, and one can use basis functions {|fc, t, x, &)} 
with multiple k for each physical orbital, which results in 
a very accurate basis. On the other hand, the MTH or¬ 
bitals completely neglect the interstitial charge and mag¬ 
netisation densities. 

The discussion on the advantages and disadvantages 
of these projection schemes can be extended to the eval¬ 
uation of exchange parameters. It is worth to mention 
that the problem is not entirely the same as that de¬ 
scribed in the previous paragraphs. For instance, for 
rare-earth metals, the 4/ electrons are exposed to pro¬ 
nounced many-body effects, which manifest themselves 
in the formation of multiplet features in the spectral func¬ 
tion. However, it is primarily the 5 d electrons which par¬ 
ticipate in the “communication” between different atoms 
by means of exchange interaction. Therefore, the local 
orbitals to use in DMFT do not have to coincide with 
those used for the Jij calculation. It was recently shown 
that the Wannier functions can be employed for the cal¬ 
culations of exchange parameters. 31 This solution is com¬ 
monly used in combination with the plane-wave-based 
DFT codes, since the basis functions are completely itin¬ 
erant. In LMTO one does not necessarily have to go via 
this step, because the LMTO basis is atomic-centered 
and sufficiently localised. The advantage of using the 
two projections schemes outlined above is that the cre¬ 
ation of the basis is less ambigous and perfectly inte¬ 
grated within the electronic structure code. When focus¬ 
ing on the exchange parameters, one can use ORT local 
orbitals without any restriction on the number of k tails. 
In fact, when one creates {|i,£,cr)} uniquely for defining 
a physical separation of the space, there is nothing that 
prevents one from using multiple basis functions of the 
same /-character. 

Since the choice of the localised states is a somewhat 
arbitrary procedure, we should analyse how it affects the 
final results. The MTH projection totally neglects the 
interstitial contribution to the magnetisation. This fact 
already implies that the Jij parameters obtained with dif¬ 
ferent projections will not be the same, simply because 
the magnetic moments are defined differently. On the 
other hand, large values of the interstitial magnetisation 
might be an indication that the system is too itinerant 
to be studied by means of Eq. (1). Hence, if the results 
strongly depend on the projection, the very meaningful¬ 
ness of the is doubtful in this case. 

To better illustrate the differences related to the local 
basis, in the upper panel of Fig. 2 we report the computed 
exchange parameters for hep Gd (upper panel) and fee 
Ni (lower panel). The computational setup of Gd is dis¬ 
cussed in section IV. The data reported here are obtained 
from DFT simulations in local spin density approxima¬ 
tion (LSDA), and for Gd we treated the 4/ electrons as 
non-hybridizing core states. The magnetic moment car¬ 
ried by the valence electrons inside the muffin-tin sphere 


is estimated for Gd to be about 0.61 /ig, in addition to 
the 7 i±b coming from the 4/ core electrons. Moreover, 
there is a 0.17 /ib per atom of the magnetisation in the 
interstitial. The latter contributes to the ORT-projected 
magnetic moment, since the underlying wavefunctions 
leak out from the MT-spheres (see Fig. 1, middle panel), 
but is absent in the MTH case. These differences among 
the local moments represent the primary reason for the 
differences between the computed J^-parameters, as seen 
in Fig. 2. Due to the truncation of the wavefunctions, 
present in MTH projection, the nearest neighbor hop¬ 
ping integrals are expected to be the most renormalised. 
The same conclusion holds for the exchange couplings, 
as confirmed by our results. Interactions with more dis¬ 
tant neighbours in the two projection schemes are not so 
different, since the detailed shape of the wavefunctions 
becomes less important. It is fundamental to note that 
in the entire range of interatomic distances, the sign and 
the relative magnitude of interactions is defined unam¬ 
biguously, i.e. does not depend on the projection scheme. 

Excluding the 4/ electrons from the consideration, 
magnetisation density in Gd is quite itinerant, which 
manifests itself in a substantial moment in the intersti¬ 
tial. Contrary to that situation, magnetisation density 
in 3d metals is much more localised around the atoms. 
Hence different projections are supposed to yield very 
similar results. For instance, in the lower panel of Fig. 2 
we show calculated Jij parameters in fee Ni. We have 
chosen rather large, almost touching radii of the MT- 
spheres in order to decrease the size of the interstitial 
region. In this case the amount of the magnetisation, 
which is not attributed to any atom is only about -0.02 
He- Magnetic moment inside the sphere (MTH projec¬ 
tion) was computed to be 0.63 he- An inspection of 
Fig. 2 reveals that different definitions of the localised 
basis sets do not significantly influence the results, even 
in the present case, where the value of a magnetic mo¬ 
ment is relatively small. 


III. EXCHANGE INTERACTIONS IN A 

WEAKLY CORRELATED METAL: BCC IRON 

The first system we choose to investigate using the 
DMFT scheme is bee iron. Its electronic structure is well- 
known and has been studied by means of DFT+DMFT 
methods by several research groups. 22,32 34 In these stud¬ 
ies it was pointed out that iron can be treated as a 
moderately correlated material. Therefore, one can re¬ 
sort to a perturbative solver, such as SPTF with val¬ 
ues of U and J parameters equal to 2.3 and 0.9 eV, 
respectively. The static part of the self-energy, aver¬ 
age over all orbitals per each spin channel, was used as 
a double-counting correction. With this choice the ex¬ 
change splitting coming from the DFT part is affected 
only through the exchange-correlation functional. The 
temperature was set to T=300 K and 1024 Matsubara 
points were sufficient to converge the self-energy E up to 
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FIG. 2. Upper panel: inter-site exchange parameters in 
hep Gd, extracted from the LSDA, where 4/ electrons were 
treated as core. Lower panel: inter-site exchange parameters 
in fee Ni, extracted from LSDA simulations. In each panel 
two sets of results were obtained by means of the two differ¬ 
ent projection schemes ORT and MTH. 


the mRy. Both MTH and ORT projections were stud¬ 
ied in DFT+DMFT and lead to similar electronic struc¬ 
tures 22 . For the present study we used two sets of basis 
functions with different k to describe the 3d states, thus 
MTH projection was employed for the solution of the ef¬ 
fective impurity problem. After that, the Jij couplings 
were computed using ORT orbitals, taking advantage of 
an extended basis set. However, as discussed above, for 
transition metals the exchange parameters obtained by 
means of the two projection schemes are very similar. 

We have extracted the exchange parameters from the 
converged LSDA, LSDA+U and LSDA+DMFT calcu¬ 
lations. The results of these calculations are shown in 
Fig. 3. Overall, the obtained Jij couplings are quite sim¬ 
ilar, at least for the LSDA and LSDA+DMFT computa¬ 
tional schemes. The LSDA+U approach deviates notice¬ 
ably, at least for the nearest-neighbor (NN) interactions, 
from the other two methods. 


The most drastic difference concerns the 2nd NN in¬ 
teraction (J 2 ), which is significantly suppressed in the 
case of LSDA+U. A more detailed analysis of the indi¬ 
vidual orbital contributions to this exchange parameter 
reveals drastically different contributions from 3d states 
with Eg and T^g symmetry, in particular we find that 
J '2 9 Eg is AFM, while J+ 29 T2a is FM. This result can 
also be explained by considering the overlap of the 3d 
wavefunctions. 35 The lobes of the T^g states of the central 
atom and its 2nd NN form a 90 degree angle, leading to 
a relatively small orbital overlap and, therefore, FM cou¬ 
pling. The Eg orbitals are pointing towards each other 
and the exchange interaction is expected to be AFM. As 
an example, for LSDA we obtain J Eg Eg =-0.06 mRy and 
J^ 2a ~ T2a = 0.62 mRy. Interestingly, in the LSDA+U ap¬ 
proximation the local exchange splitting of the E g states, 
being about 0.21 Ry, is enhanced with respect to the 
LSDA value of 0.16 Ry. The opposite tendency is found 
for the T^g orbitals (0.06 Ry versus 0.13 Ry, respectively). 
Since the local exchange splitting A enters Eq. (2) twice, 
any changes of this quantity strongly influence the effec¬ 
tive exchange parameter. As a result, the J Eg Eg com¬ 
pletely compensates the J E2a T2a contribution in the case 
of LSDA+U. We checked that further increase of U leads 
to the overall change of the J 2 sign, because the AFM 
contribution starts to prevail. 

The data shown in Fig. 3 was produced using the same 
U and J parameters for LSDA+U and LSDA+DMFT. 
The latter approach, however, contains also the 3d 
screening on the effective on-site potential, which is to¬ 
tally neglected in the static mean-field LSDA+U. There¬ 
fore, for the same strength of U one can expect that 
the results of LSDA+DMFT will be somehow closer to 
that of the conventional LSDA calculation. This expec¬ 
tation is in line with the results of our numerical simu¬ 
lations. We can also use the computed exchange param¬ 
eters to evaluate the ordering temperature in the mean- 
field approximation (MFA): T^ FA = Jo, where Jo = 
Joj- We obtain a value of 925 K in LSDA, 905 K 
in LSDA+U and 840 K in LSDA+DMFT. The experi¬ 
mental value of T c is about 1045 K, but any agreement 
in absolute numbers should be regarded as somewhat 
fortunate. In fact we should mention that, in general, 
the MFA gives higher values of T c when compared with 
Monte-Carlo simulations. This is mainly due to that nu¬ 
merical estimations are based on the assumption that 
there is no short-ranged order above the critical tem¬ 
perature and there is no temperature dependence of the 
exchange parameters. However, both approximations are 
in general rather questionable. In fact recent work sug¬ 
gests that the inter-atomic exchange does depend on the 
temperature driven magnetic configuration. 36 

The differences in the interactions for the more distant 
neighbours, shown in Fig. 3, are difficult to distinguish. 
It is known that RKKY interactions 37 are expected to 
have a R~ 3 asymptotics, where Rij is the interatomic 
distance. Hence, in order to illustrate the long-ranged 
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FIG. 3. Exchange parameters in bee iron (in mRy) as a func¬ 
tion of the interatomic distance (in the units of the lattice 
constant). Positive sign corresponds to the ferromagnetic cou¬ 
pling. ”ORT” projection scheme was employed in the present 
calculation of the inter-atomic exchange. DMFT results were 
produced using the SPTF solver. 


character of the magnetic couplings, Fig. 4 contains the 
exchange parameters with the first 99 shells, multiplied 
with Rfj . These results were obtained with a very dense 
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FIG. 4. Exchange parameters Jij in bcc iron, factorised by 
inter-atomic distance R%. 


k-point grid, consisting of 96 3 points in the entire Bril- 
louin zone (BZ). This is an unnecessary large number 
for the accuracy of the local quantities, but is, however, 
necessary for that of the inter-site Green’s function. 5 

Inspection of Fig. 4 reveals that indeed LSD A and 
LSDA+DMFT results are extremely similar. It is worth 
mentioning that this is also related to the fact that we 
use SPTF solver for the DMFT problem, which intrinsi¬ 
cally forces moderate correlation effects, but this is what 
is expected for bulk bcc Fe. One can see that the effect 
of dynamical correlations is such that the Jij parameters 


are renormalised with respect to the LSDA values. The 
renormalisation is, however, anisotropic and is dependent 
on the direction of the interatomic bond. Regarding the 
comparison of the absolute values, one also has to bear 
in mind that the calculated magnetic moment per Fe is 
about 2.2 iib in LSDA, whereas it is 1.95 /r b and 2.15 
Hb in LSDA+U and LSDA+DMFT, respectively. In the 
considered form of the Hamiltonian (Eq. (1)) it means 
that this difference is effectively contained within the Jij 
parameters. 

In order to quantify the differences in the obtained ex¬ 
change parameters, we have simulated the frozen magnon 
spectra, expressed through their Fourier transform J(q) 
(see e.g. Ref. 5). Fig. 5 shows the simulated magnon 
dispersions along [100] direction together with the in¬ 
elastic neutron scattering data. For small | q | the results 
obtained with different electronic structure methods are 
quite similar. Most of the differences become evident at 
higher momenta, but these regions are difficult to access 
experimentally. Within the considered range of g-space, 
the agreement with experimental data is quite good. 



FIG. 5. Magnon dispersion along [100] direction in the BZ 
of bcc Fe, simulated using the computed Jij -parameters. All 
calculations were performed for T=300 K. Experiment: room 
temperature inelastic neutron scattering data (from Ref. 38) 
are shown with circles. 


Recently Mazurenko et al. have shown that in the 
Fermi-liquid regime the total exchange interaction Jq for 
a given site is supposed to be rescaled with an effective 
mass enhancement factor Z. 39 In general, our results are 
in line with their conclusion, since we also report an over¬ 
all decrease in the magnitude of the magnetic couplings 
once correlations are introduced. On the other hand, 
within the LSDA+DMFT approach the self-energy is dif¬ 
ferent for both spin channels and so is the Z factor. More¬ 
over, during the charge self-consistency the unperturbed 
DFT Hamiltonian gets also modified. Most importantly, 
our analysis of the next NN coupling in bcc Fe demon¬ 
strated that if there are several exchange paths, the mod- 
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ifications of the exchange parameters can be highly non¬ 
trivial. Thus the relation between Jo obtained in LSDA 
and that in LSDA+DMFT is more complex. 


IV. STRONGLY CORRELATED SYSTEM: 

NICKEL OXIDE 

Next, we have applied our implementation to the 
canonical example of transition metal oxides (TMOs), 
namely NiO. Magnetic interactions in these materials 
are rather short-ranged, since they are insulating, hence 
considering a few neighbouring shells is usually suffi¬ 
cient to describe their spin-wave spectra. 40,41 . In the 
literature, the exchange parameters have been success¬ 
fully described by means of the LSDA+U approach. 42 44 
Concerning the details of the electronic structure the 
LSDA+U is only moderately successful. Magnetic mo¬ 
ments and band gap are well described in these materials, 
but the comparison between experimental photo-electron 
spectra and theoretical one-particle excitation spectra is 
far from satisfactory 45 . Hence, in order to reproduce 
both valence band 21 and magnon 17 spectra within the 
same framework one must use the more sophisticated 
LDA+DMFT method, as we argue here. 

Since the correlation effects are known to be very pro¬ 
nounced for these compounds, we could not use the SPTF 
solver for the impurity problem, but we preferred the ED 
solver. Differently from Ref. 17 we performed DMFT 
simulations starting from the LSDA solution. In this 
way the spin-polarisation already exists at the level of the 
DFT part of the problem, i.e. before the solution of an 
effective impurity problem. We employed full-localised 
limit DC, which is appropriate for the systems close to 
atomic localisation. 46 At convergence, the exchange split¬ 
ting is due to both the effects of the exchange-correlation 
functional and those of the dynamical self-energy. 47 

The DFT+DMFT simulations with ED were obtained 
with the same parameters used in Ref. 21. The fitting of 
the hybridization function was done with two bath sites 
per each 3d orbital. The V and J parameters for Ni 3d 
states were set to 8.2 and 0.95 eV, respectively. In order 
to reach the saturation of the local magnetisation, the 
temperature was lowered down to 80 K, thus leading to 
an increase of the number of Matsubara frequencies up 
to 6144 points. In the following we only present results 
obtained with the ORT projection. 

According to the results of LSDA+DMFT calculation, 
NiO is characterised by a band gap of about 3.5 eV and 
1.85 ii b of spin moment per Ni ion. The LSDA+U sim¬ 
ulations gave similar results. The calculated magnetic 
moments are in fair agreement with the experimental 
value of 1.90 /i#. 48 The calculated exchange parameters 
in NiO, obtained by means of different techniques are 
listed in Table I. 

It is known that the LSDA gives a rather poor de¬ 
scription of the electronic structure of the TMOs. In 
the present calculation we obtained a tiny band gap of 


TABLE I. Magnetic interactions in NiO, computed by means 
of various ab initio techniques. Exchange parameters Jn 
(in mRy) denotes the interaction with the iV-th coordination 
shell. DMFT results were obtained using ED as an impurity 
solver. 


Computational setup 

Ji 

j 2 

LSDA 

0.04 

-1.58 

LSDA+DMFT 

-0.003 

-0.48 

LSDA+U 

-0.002 

-0.50 

LSDA+U (U=8 eV) [Ref. 42] 

0.004 / 0.0 

-0.53 

Exp. # I [Ref. 41] 

-0.051 

-0.637 

Exp.#2 [Ref. 49] 

0.051 

-0.67 


about 0.5 eV, much smaller than the experimental esti¬ 
mate. Since the system is an insulator, at least the un¬ 
derlying exchange mechanism is expected to be correct, 
i.e. super-exchange is expected. Note that the situation 
is much worse for other TMOs, which are predicted to 
be metals and hence electrons at the Fermi surface par¬ 
ticipate in the magnetic coupling, resulting in a double 
exchange mechanism. As a result, the sign of the leading 
Jij interaction in NiO is correct, but its absolute value is 
way too overestimated in LSDA, which is a well-known 
fact. 

Overall, the J^-parameters obtained in DMFT are 
very similar to those obtained in LSDA+U, both in the 
present and prior calculations. 42 ” 44 The magnetic inter¬ 
actions in this class of systems are defined primarily by 
the virtual electron hoppings and relative positions of the 
electron levels. Both of these quantities are described rel¬ 
atively well already on the level of LSDA+U. NiO is a 
particularly illustrative example, because there is a single 
d orbital, which participates in the dominant J 2 coupling. 
The latter, according to the super-exchange mechanism, 
is expected to be proportional to t 2 /U, where t is an ef¬ 
fective inter-site hopping integral. 50 To test the relevance 
of the super-exchange mechanism in a quantitative way, 
we show in Fig. 6 the calculated J 2 coupling as a function 
of the V parameter. 

An inverse proportionality between J 2 and U is recog¬ 
nisable. If the exchange is given exactly by the term t 2 /U 
one would expect that + 2 * 1 / should result in a constant 
value, t 2 . In order to test this we show (—J 2 * U) in 
the inset of Fig. 6. One may note that the so obtained 
value of t 2 is essentially constant and more or less inde¬ 
pendent on U, illustrating that the theory captures cor¬ 
rectly the super-exchange mechanism of NiO. The results 
of Fig. 6 also imply that the effective hopping parame¬ 
ters extracted from LSDA+U and LSDA+DMFT agree 
well with each other. On the other hand, a weak but 
non-negligible dependence of the effective hopping with 
respect to U is discernible. The reason for that is that the 
calculations are done in the charge self-consistent man¬ 
ner, which implies the change of the DFT potential and, 
therefore, t values will depend on the U value used in the 
calculations. In fact, in these calculations the V value 
was varied in an unrealistically wide range. Considering 
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V. RARE-EARTH METAL: HCP Gd 



FIG. 6. Next NN exchange parameter in NiO for vari¬ 
ous values of Hubbard U calculated within LSDA+U and 
LSDA+DMFT. Inset: The same plot for (—J 2 * U). 


a narrower window of more reasonable values, one can 
approximate t to be essentially constant. 

Hence we can conclude that the improvements of the 
use of DMFT are quite moderate for this system, mean¬ 
ing that static local correlations capture most of the es¬ 
sential modifications from the LSD A treatment. This is 
to some extent surprising since the features of the elec¬ 
tronic structure of the DMFT and LSDA+U approxima¬ 
tions are very different. However, one may argue that 
in order to get the correct value of super-exchange one 
needs a correct value of the hopping parameter and the 
gap, and it seems LSDA+U does a good job in providing 
these properties accurately. We suspect that the present 
conclusion holds for most wide gap insulators with large 
localised magnetic moments. 

We also report a fair agreement with experimental 
data, but a few words about these data are needed. From 
Table I we see that different experiments lead to slightly 
different values of the Jij -parameters, and there is not 
even consensus about the sign of the NN coupling. The 
experimental data reported in Table I are obtained by fit¬ 
ting the measured magnon spectra to an effective Heisen¬ 
berg model. The precision of this procedure is question¬ 
able, especially if one consider that Refs. 41 and 49 report 
on the studies that are forty years old. From the numer¬ 
ical side, as seen from Fig. 6, by choosing a bit smaller 
V (about 6.8 eV) one can obtain J 2 = 0.62 mRy. With 
this respect, we consider our theoretical values in good 
agreement with experimental data. Our analysis leads 
to the AFM NN coupling, but the absolute value of this 
interaction is very small. 


The 4/-shell of heavy rare-earth metals is sufficiently 
localised and essentially does not participate in the chem¬ 
ical bonding. However, if treated at the LDA level, these 
states acquire finite dispersion and are placed at the 
Fermi level in contradiction with existing experimental 
data. In particular, in hep Gd the LSDA approxima¬ 
tion leads to an appearance of antiferromagnetic 4/ — 5 d 
exchange interaction, which destabilises the true ferro¬ 
magnetic ground state (see e.g. Ref. 51). One possibility 
to overcome this issue is to use the LSDA+U approach, 
which results in a static U-correction, shifting the 4/ 
states away from the Fermi level. 46 It was also shown 
that treatment of 4/ states as purely non-hybridising 
core levels, improves the description of inter-site mag¬ 
netic couplings. 52,53 There was, however, a discrepancy 
in the calculated J^-parameters, which can be explained 
by the use of different computational codes. Here we 
have investigated the electronic and magnetic structure 
of metallic Gd by describing 4/ states within HIA. The U 
and J parameters were respectively set to 8.2 eV and 0.7 
eV, in agreement with prior studies. 46,54 Given that our 
localised basis is different from those used in these other 
studies, a one-to-one correspondence of the exchange pa¬ 
rameters can not be expected. 

The HIA can be seen as a trivial solution in DMFT, 
where one performs an ED calculation by removing all 
the hybridization of the local orbitals with all the other 
states. However, one is still able to foster multiplet peaks, 
which are measured experimentally. The calculated elec¬ 
tronic structure from the LSDA+DMFT calculations is 
depicted in Fig. 7, where it is compared to experimental 
data. 

An inspection of the plot reveals that the 4/ states 
form two distinct groups of multiplets, located around 
—8 eV and 4 eV relative to the Fermi level. These two 
groups of states have purely opposite spin polarisations. 
The agreement with experimental data is acceptable, but 
can be improved by a slightly different choice of U and J, 
which, however, will not affect the overall picture of the 
magnetic properties of the system. Moreover, the photoe¬ 
mission measurements were done at room temperature, 
above T c , where the electron levels are spin-degenerate. 
We note here that because the ground state of Gd has 
a filled spin-up shell and empty spin-down shell, the ex¬ 
citations of this configuration become particularly sim¬ 
ple and actually quite similar to results from LSDA+U. 
However, for essentially all other rare-earths, either as 
compounds or in elemental form, the LSDA+U method 
does not reproduce the measured multiplet structure. 

Regarding the exchange interaction parameters, the 
values are shown in Fig. 8. We found that different treat¬ 
ment of 4/ states (by means of HIA or simply treating 
them as part of the core) does not give any qualitative 
changes. The most dramatic change is associated with 
the NN couplings. There is an anisotropy of the in-plane 
and out-of-plane exchange parameters, which is a result 
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FIG. 7. Total and /-projected spectral function in hep Gd 
as computed in LSDA+DMFT method within the Hubbard-I 
approximation. The Fermi level is at zero energy. X-ray pho¬ 
toelectron spectroscopy (XPS) and bremsstrahlung isochro- 
mat spectroscopy (BIS) data is adapted from Ref. 55. 
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FIG. 8. Inter-site exchange parameters in hep Gd, extracted 
from the ground states with different treatment of 4/ elec¬ 
trons. The MTH-projection scheme was adopted to produce 
both sets of results. 


of an experimentally measured reduced c/a ratio (equal 
to 1.59 in the current setup). This difference is enhanced, 
when the 4/ electrons enter explicitly in the calculation, 
which seems to be related to a more non-uniform dis¬ 
tribution of the electron density. The rest of the Jij pa¬ 
rameters are quite similar overall: their signs and relative 


magnitudes are identical. The spin magnetisation, asso¬ 
ciated with the most magnetically “active” 5 d electrons, 
is about 0.48 hb- This results in an induced valence 
band moment of 0.68 Hb /&tom (including the interstitial 
and 6s/6p MT contributions). The latter is almost inde¬ 
pendent on the 4/ treatment employed. The computed 
exchange parameters are of long-ranged nature, clearly 
showing RKKY oscillations. Due to this fact, we em¬ 
phasize that in order to have a converged estimate of T c 
or magnon spectrum (not shown here), we have to em¬ 
ploy interactions with at least 23 first neighbouring shells. 
Our results are in good agreement with prior studies, 
given that the computational schemes are different (e.g. 
an adoption of the atomic sphere approximation). 52,53 
Most importantly, our results are obtained with a full- 
potential electronic structure code, that can treat open 
structures, surfaces, etc. 


VI. CONCLUSIONS 

We report here on the implementation of interatomic 
exchange parameters for correlated electronic structures 
and analyse the exchange mechanisms of three classes of 
correlated electron systems; the weakly correlated tran¬ 
sition metal elements, the strongly correlated transition 
metal oxides and the localized electron systems repre¬ 
sented by the rare-earth element Gd. We find that de¬ 
tails in how electron correlations are included influence 
the exchange interactions in different amounts, depend¬ 
ing on which class one considers. To be more specific, we 
analyse the modifications of exchange couplings in weakly 
correlated systems, like bee Fe, as compared to the re¬ 
sults of conventional LSDA. Since the correlation effects 
act differently on each orbital, their contributions to the 
total exchange coupling can be substantially different. In 
the case of transition metal oxides, static on-site corre¬ 
lations are shown to be the most important ones for an 
adequate description of the exchange couplings in these 
systems. Finally we have applied the method to hep Gd, 
where correlated orbitals are of / nature, but more de¬ 
localised spd states actually participate in the inter-site 
magnetic interactions. We show that the J^-parameters 
are barely affected by different treatment of 4/ states, 
because the corresponding states hybridise weakly with 
the rest. 

Regarding the comparison with available experimental 
data, it is clear that the LSDA offers an adequate de¬ 
scription of the magnetic properties of weakly correlated 
metals, whereas the LSDA+U is suitable for wide gap 
insulators. However, only the DMFT-based calculations 
provided us a systematically good agreement for all the 
materials studied. This fact underlines that the inclu¬ 
sion of electron correlations beyond that of LSDA is in 
general important. 

Using the developed machinery one can study the tem¬ 
perature dependence of the magnetic excitations as one 
approaches T c . A few words have to be said about the 
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calculation of the T c . In DMFT one does not take into 
account all possible finite temperature effects. Even by 
using the most accurate quantum Monte-Carlo (QMC) 
solver with an approximate matrix of Coulomb inter¬ 
actions (density-density approximation), there are only 
longitudinal fluctuations, which contribute to the desta¬ 
bilisation of the long range magnetic order (and hence 
define the T c ). This is opposite with respect to the linear 
response technique, which we adopt, since here one deals 
with the transverse spin fluctuations solely. Therefore, 
the estimates of the T c obtained with s as calculated 
here, in combination with MFA or Monte Carlo simu¬ 
lations, are not expected to be the same as the results 
obtained using only DMFT and QMC. 

Recently an extension of the method 3,12 was proposed, 
which permits to extract magnetic couplings from non- 
collinear spin configurations. 36 This is a promising ap¬ 
proach, especially if one wants to model the magnetic 
properties at finite temperature, where the partial reori¬ 
entation of the spins is present. Moreover, the employed 
linear-response approach (at least in the current form) as¬ 
sumes that the dynamics of the spin waves is adiabatic, 
i.e. the unperturbed electronic structure is used to ex¬ 
tract effective exchange parameters. In general case, the 
adiabaticity is not guaranteed and magnetic excitations 
should be regarded as truly out-of-equillibrium processes. 
The formulae for calculating the magnon spectra in such 
case have been recently proposed by Secchi et al , 56 How¬ 


ever, for the moment the method can only be applied 
to the model systems and its reconciliation with the ab 
initio techniques is far from being at hand. 

Many magnetic systems exhibit interesting proper¬ 
ties due to the presence of other types of magnetic 
interactions, originating from the spin-orbit coupling. 
Among those there are Dzyaloshinskii-Moriya (DM) 
interaction 57,58 and magnetocrystalline anisotropy (MA). 
An interplay between different terms gives rise to 
non-collinear spin structures, observed in multiferroic 
materials, 59 magnetic ions on surfaces, 60 and some low- 
symmetric TM alloys. 61 Both DM and MA can be cal¬ 
culated using the linear-response-like techniques. 10,62-65 
However, to our knowledge, the effects of dynamical cor¬ 
relations on the anisotropic parameters were never stud¬ 
ied so far. Our implementation provides an ideal base for 
conducting an investigation of this problem. 
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